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ABSTRACT 

We use a simulation performed within the Constrained Local UniversE Simulation 
(CLUES) project to study a realistic Local Group-like object. We employ this group 
as a numerical laboratory for studying the evolution of the population of its subhaloes 
from the point of view of the effects it may have on the origin of different types of 
dwarf galaxies. We focus on the processes of tidal stripping of the satellites, their 
interaction, merging and grouping before infall. The tidal stripping manifests itself 
in the transition between the phase of mass accretion and mass loss seen in most 
subhaloes, which occurs at the moment of infall on to the host halo, and the change 
of the shape of their mass function with redshift. Although the satellites often form 
groups, they are loosely bound within them and do not interact with each other. The 
infall of a large group could however explain the observed peculiar distribution of the 
Local Group satellites, but only if it occurred recently. Mergers between prospective 
subhaloes are significant only during an early stage of evolution, i.e. more than 7 Gyr 
ago, when they are still outside the host haloes. Such events could thus contribute to 
the formation of more distant early type Milky Way companions. Once the subhaloes 
enter the host halo the mergers become very rare. 

Key words: methods: iV-body simulations - galaxies: Local Group - galaxies: dwarf 
- galaxies: fundamental parameters - galaxies: kinematics and dynamics - cosmology: 
dark matter 
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1 INTRODUCTION 

According to the favoured ACDM model galaxies are formed 
in a hierarchical way (White & Rees 1978). Small galaxies 
form first and then merge into larger structures. In this sce- 
nario the surviving dwarf galaxies can be viewed as tracers 
of the early Universe. Understanding their properties and 
evolution is crucial to understanding the ACDM Universe 
itself. 

The Local Group (LG) is the closest and best stud- 
ied object of the extragalactic scale (see van den Bergh 
1999 for a review). Apart from the two large spirals, the 
Milky Way (MW) and Andromeda (M31), it is populated 
by several tens of dwarf galaxies which by their morpholog- 
ical properties can be divided into irregulars, ellipticals and 
spheroidals (Mateo 1998). They are known as the classical 
dwarfs. Recently we have witnessed many new discoveries 
of ultra faint satellite galaxies in the halo of the MW (e.g. 



Sakamoto & Hasegawa 2006; Zucker et al. 2006; Belokurov 
ct al. 2008). The dynamical properties of these objects are 
still very poorly studied, but their alleged extremely high 
mass to light ratios (Simon & Geha 2007) suggest that they 
might be a rather different class of objects. 

One could expect that an early Universe progenitor 
of a dwarf galaxy should be similar for all types of the 
present classical dwarfs (Kravtsov, Gnedin & Klypin 2004) 
and would consist of a dark matter halo of mass around 10 9 
Mq. This amount of mass would allow it to accrete gas even 
in the highly ionized environment of the early Universe. The 
conservation of angular momentum would require that the 
gas and the stellar component formed a disk rather than a 
spheroid. Such an object would then undergo evolution due 
to baryonic processes such as cooling, star formation, super- 
nova feedback, photoevaporation of the gas and heating by 
the cosmic ultraviolet background. It remains to be seen if 
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Figure 1. Overview of the simulated LG. M31 is located in the 
centre of the image, MW at the top right. The image size is 3 
h~ 1 Mpc. 

these or rather the environmental, purely gravitational pro- 
cesses were the dominant factor in shaping dwarf galaxies. 
In this work we focus on the latter which can be modelled 
with simulations following only the dark matter component. 

The studies of stellar populations and kinematics of 
dwarf galaxies provide some general hints (see Tolstoy, Hill 
& Tosi 2009 for a recent review) concerning the possible 
environmental effects. The late type dwarfs (dwarf irregu- 
lars) possess rich star formation histories, high amounts of 
hydrogen in stellar disks and significant angular momenta. 
They are usually isolated objects (except for the Magellanic 
Clouds). These facts suggest that they evolved in isolation 
only by the internal baryonic processes. On the other hand, 
the early type dwarfs usually lack gas and do not show any 
signs of recent star formation. They also possess very low 
or no angular momentum. They typically are close compan- 
ions to the MW and M31, but even those isolated ones like 
Tucana or Cetus could have evolved in the vicinity of the 
large spirals and been ejected recently to more distant orbits 
(Gill, Knebe & Gibson 2005; Sales et al. 2007a; Ludlow et 
al. 2009). They have high mass-to-light ratios from about 
ten to a few hundred solar units. They could have evolved 
from a disk like those present in dwarf irregulars by strong 
gravitational interactions. These interactions could take the 
form of mergers or tidal forces. 

The latter idea has been developed into the so-called 
tidal stirring scenario (Mayer et al. 2001). Using N-body 
simulations it has been demonstrated that indeed a trans- 
formation from a disk to the dwarf spheroidal is possible 
by pure tidal interactions with a host galaxy. Mayer et al. 
(2001) have shown that low surface brightness disks would 
produce dwarf spheroidal galaxies we observe today, while 
high surface brightness disks would lead to the formation of 
dwarf ellipticals. Depending on the adopted star formation 
history this scenario is able to reproduce dwarf spheroidals 
with moderate mass-to-light ratios (Klimentowski et al. 
2007, 2009) as well as strongly dark matter dominated ones 
if the gas is expelled early on (Mayer et al. 2007). 

Mergers and interactions between dwarfs could be an- 



Table 1. The properties of the LG: distance and relative velocity 
between the two main members, the velocity dispersion around 
the Hubble flow, the overdensity within a sphere of 7 Mpc centred 
on the LG, distance from the LG to the Virgo cluster and the 
virgocentric flow. The second column gives the values from the 
simulation. For comparison, the third column provides also the 
observed values. The references for the observed values in the last 
column arc: 1 - Karachcntscv et al. (2004); 2 - van der Marel & 
Guhathakurta (2008); 3 - Tikhonov & Klypin (2009); 4 - Hudson 
(1993); 5 - Fouque et al. (2001); 6 - Tonry et al. (2000). 



Property 


Simulated 


Observed 


Ref. 


d (M31-MW) 


0.91 h' 1 Mpc 


0.77 Mpc 


1 


V Icl (M31-MW) 


-193 km s" 1 


-130 km s" 1 


2 


<r H (r < 7 Mpc) 


100 km s _1 


90 km s" 1 


3 


Sp/p (r < 7 Mpc) 


0.82 


0.8 


4 


d (LG- Virgo) 


10.9 h^Mpc 


18.0 Mpc 


5 


Virgocentric flow 


-255 km s" 1 


-409 km s" 1 


6 



other channel for the formation of early type dwarfs since 
they could lead to a very strong evolution (e.g. Knebe et al. 
2006; Angulo et al. 2009), as in the case of larger galaxies 
(e.g. Springel, Di Matteo & Hernquist 2005). It has also 
recently been proposed that dwarf galaxies might be ac- 
creted on to their hosts in groups which may explain e.g. 
the particular distribution of dwarfs around the MW and 
M31 (Libeskind et al. 2005; Metz, Kroupa & Jerjen 2007, 
2009a; D'Onghia & Lake 2008; Li & Helmi 2008; Metz et al. 
2009b). In this work we discuss this second scenario of the 
formation of early type dwarf galaxies using a cosmological 
JV-body simulation of a LG. Our purpose is to study merg- 
ers and interactions of haloes which end up as subhaloes. We 
also consider the infall of these objects together as a group 
of small haloes. 

The paper is organized as follows. Section 2 contains 
the description of the simulation used in this analysis and 
the halo finding algorithm which provides the basis for this 
study. In section 3 we characterize the main properties of 
the subhalo population of the two largest haloes; we discuss 
their mass functions, their survival times and the evolution 
of their masses. In section 4 we study the behaviour of sub- 
haloes in groups; we provide the statistics of groups around 
the two largest haloes, describe the mass functions of the 
largest groups and follow their history. The effect of mergers 
and interactions between subhaloes is discussed in section 5. 
The discussion follows in section 6. 



2 THE SIMULATION 

We analyze a constrained dark matter simulation of the LG. 
This simulation is part of the CLUES projedQ, a collabo- 
ration whose main goal is to produce realistic cosmological 
simulations of the Local Universe by imposing observational 
constrains on the mass and velocity fields of the initial ran- 
dom Gaussian fluctuation realizations. For this simulation 
we have used a box of 64h~ Mpc size assuming a spatially 
flat cosmological model with WMAP3 parameters (Spergel 

1 http://clues-project.org 
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Figure 2. The cumulative mass function of subhaloes found within one virial radius of two most massive haloes normalized by the 
volume of that region. The left panel shows N(> M) at redshift z = 0. The right panel presents the mass function of the same haloes 
traced back to redshift 2 = 1. The upper axis indicates the number of particles in a halo. Vertical dotted lines indicate our completeness 
limit of 100 particles. 



Table 2. The properties of the five most massive haloes identified in the simulation. The columns list the assigned name of a halo, its 
virial mass, the virial radius, the number of particles inside the virial radius and the spin parameter (according to Peebles definition). 



Name M vir [lO 11 fe" 1 M ] rvir^^kpc] N A 



Andromeda (M31) 


5.69 


173.5 


2239732 


0.0686 


Milky Way (MW) 


4.62 


162.0 


1821419 


0.0607 


3rd halo 


2.66 


134.7 


1045831 


0.0722 


4th halo 


2.24 


127.2 


883066 


0.0219 


5th halo 


2.01 


122.7 


738875 


0.0420 



et al. 2007): Q, m = 0.24, fi b = 0.042, fi A = 0.76, the Hub- 
ble constant h — 0.73, the normalization as = 0.75 and the 
slope n = 0.95 of the power spectrum. 

First, a constrained density field on a grid of 256 3 mesh 
points was obtained applying the Hoffman & Ribak (1991) 
algorithm for generating constrained realizations of Gaus- 
sian random fields. As observational constraints we have 
used the radial velocities of galaxies drawn from the MARK 
III catalogue (Willick et al. 1997), Surface Brightness Fluc- 
tuation survey (Tonry et al. 2001) and the local volume 
galaxy catalog (Karachentsev et al. 2004) as well as the posi- 
tions of nearby X-ray selected clusters of galaxies (Reiprich 
& Bohringer 2002). The algorithm has been described in de- 
tail in Zaroubi, Hoffman & Dekel (1999), Kravtsov, Klypin 
& Hoffman (2002) and Klypin et al. (2003). With this algo- 
rithm to calculate the initial conditions the resulting sim- 
ulation contains the main features which characterize the 
Local Universe. In the large simulation box the Virgo as 
well as the Coma cluster and the Great Attractor are ap- 
proximately at the right positions whereas the small scale 
structure is essentially random. A smaller box like the one 
we discuss here contains an object which can be identified 
as the Virgo cluster. 

Within our simulation box we represent the linear power 
spectrum at redshift z = 100 by N max = 4096 3 particles of 
mass ?tidm = 2.5 x 10 s h~ 1 M@. We then Fourier transform 



the constrained density field and substitute the overlapping 
Fourier modes in our otherwise random realization. At first 
we degrade the mass resolution to 256 3 particles and identify 
the position of the simulated LG at z = 0. To this end 
we start at the position of the simulated Virgo cluster and 
search for a Local Group like object at the right position. 
After identifying such an object we find all the particles 
within a sphere of radius 2h~ Mpc centred on the simulated 
object and determine the Lagrangian coordinates of these 
particles in the initial conditions. 

In a next step we resimulate the evolution of the LG 
using the full resolution (4096 3 ) within this sphere of radius 
2/i _1 Mpc. Here we follow the algorithm described in Klypin 
et al. (2001) and degrade the mass and force resolution in 
those areas that are far away from the Lagrangian region 
from which the LG forms. To this end we put concentric 
regions around the high-resolution area, each of them pop- 
ulated with particles 8 times more massive. In the end we 
have five different mass refinements ranging from 4096 3 at 
the high resolution area to 256 3 in the outer parts of the 
simulation box. Thus we simulate with very high resolution 
the evolution of the LG in the right environment. The evo- 
lution of the same region has been also simulated including 
gas dynamics and star formation (Libeskind et al. 2009). 

The simulation has been performed using the TreePM 
parallel AT-body code GADGET2 (Springel 2005). For the 



© 0000 RAS, MNRAS 000, 000-000 



4 J. Klimentowski et al. 



high mass resolution particles we used a fixed comoving 
Plummer equivalent softening of 500 h~ 1 pc at early red- 
shift and changed to 100 hT x pc physical since z = 4. For the 
rest of the mass refinement levels, we increase the Plummer 
softening by a factor of approximately two at every level. To 
follow the evolution of the LG we have stored in total 134 
outputs equally spaced in time, which translates into a time 
difference of 0.1 Gyr between consecutive snapshots. The 
overview of the simulated LG in the final output is shown in 
Fig- [Hand its properties are listed in Table[T] The properties 
of its five most massive haloes are given in Table [2] 

In order to identify haloes and subhaloes in our simu- 
lation we have run the MPI+OpenMP hybrid halo finder 
AHF (AMIGA halo finder, to be downloaded freely from 
http://www.popia.ft.uam.es/AMIGA). AHF is an improved 
version of the MHF halo finder (Gill, Knebe & Gibson 2004a), 
which locates local overdensities in an adaptively smoothed 
density field as prospective halo centres. The local poten- 
tial minima are computed for each of these density peaks 
and the gravitationally bound particles are determined. We 
stress that our halo finding algorithm automatically iden- 
tifies haloes, subhaloes, subsubhaloes, etc. For more details 
on the mode of operation and the actual functionality we 
refer the reader to Knollmann & Knebe (2009) where the 
AHF halo finder is described in detail. 

For each halo, we compute the virial radius r v j r at which 
the density M(< r)/(47rr 3 /3) drops below A v j r pb where 
pb is the cosmological background density. The threshold 
Avir is computed using the spherical top-hat collapse model 
and is a function of both cosmological model and time (e.g. 
Lokas & Hoffman 2001). For the cosmology we are using, 
Avir = 355 at z = 0. Subhaloes are defined as haloes which 
lie within the virial region of a more massive halo, the so- 
called host halo. As subhaloes are embedded within the 
density of their respective host halo their very own den- 
sity profile usually shows a characteristic upturn at a radius 
r t < r v i r , where r v i r would be their actual (virial) radius 
if they were found in isolation^ We use this 'truncation ra- 
dius' r t as the outer edge of the subhalo and hence (sub)halo 
properties (i.e. mass, density profile, velocity dispersion, ro- 
tation curve) are calculated using the gravitationally bound 
particles inside either the virial radius r v i r for a host halo or 
the truncation radius r\ for a subhalo. 

Once the halo finding is completed all haloes are traced 
back in time. For this purpose the halo must be linked to 
its progenitor in the previous simulation output. We do 
that using the following prescription. A halo progenitor is 
identified in the previous snapshot by maximizing the ratio 
Cf /(NiNj), where d is the number of common particles 
shared between the i-th halo of the current snapshot and 
the j'-th halo of the previous snapshot while Ni and Nj are 
the total numbers of particles in these haloes. 

This simple formula works surprisingly well in identify- 
ing the correct halo in the past. Unfortunately, sometimes 
the correct halo is missing in the previous snapshot because 



2 Note that the actual density profile of subhaloes after the re- 
moval of the host's background drops faster than for isolated 
haloes (e.g. Kazantzidis et al. 2004a); only when measured within 
the background still present we will find the characteristic upturn 
used here to define the truncation radius rt. 
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Figure 3. The ratio between the final mass of a subhalo at the 
end of the simulation and its mass at the infall time, as a function 
of infall time for M31 (upper panel) and MW (lower panel). 



it was not found by the halo finder. This happens, for ex- 
ample, when a smaller halo passes close to the centre of a 
larger one and due to its low density contrast the halo finder 
may not be able to identify it. To take care of this problem 
(wrong identification would lead to something which looks 
like a halo splitting in two) we also check that the mass 
of the progenitor is close to the mass of the descendant: if 
the mass ratio is smaller than 0.8 we look for the correct 
progenitor halo two (instead of one) snapshots earlier. This 
parameter was tuned to the actual properties of the halo 
finder. This procedure is applied recursively until we find a 
credible progenitor in one of the higher redshift snapshots. 
In practice, however, we never have to skip more than one 
or two snapshots. 
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Figure 4. Example of the evolution of a satellite subhalo. The solid line shows the mass of the subhalo and the dashed one its distance 
from the host galaxy divided by its evolving virial radius, as a function of time. The dotted line shows the distance of one virial radius 
from the host, which we adopt as a measure of the orbit entry point for the subhalo. During the first stage of the evolution the subhalo 
gains mass. After entering the orbit around its host it starts losing mass due to tidal stripping. 



3 PROPERTIES OF THE SATELLITES 

The first two rows of Table [5] list the properties of the two 
most massive haloes identified in the last simulation output 
corresponding to redshift z = 0. We will call the largest one 
Andromeda (M31) and the second largest the Milky Way 
(MW). In this section we study in detail the properties of 
the satellites of these largest haloes. The left panel of Fig. [2] 
shows the mass function N(> M) of all subhaloes belong- 
ing to the two most massive haloes at redshift z = 0. The 
subhaloes were identified within one virial radius of their 
host. The relation was scaled by the volume of each region 
in which subhaloes were selected. In this case the mass func- 
tions are very similar for both haloes. Note that the shapes of 
the mass functions flatten towards smaller masses signifying 
loss of completeness. This happens at masses corresponding 
to about 100 bound particles or 2.5 x 10 7 Ti" 1 M©. Unless 
otherwise stated, from now on we will consider only haloes 
with masses above this value. 

The right panel of Fig. [2] shows the cumulative mass 
function of all subhaloes found at redshift z — traced back 
to z = 1. Comparing the panels we see that with decreasing 
redshift the number density of small mass haloes increases 
and the whole slope of the relation steepens. This suggests 
that subhaloes which were originally more massive lose mass 
during their evolution. This process is due to the tidal in- 
teractions of subhaloes with their host halo and is common 
for satellites of all large galaxies (cf. Mayer et al. 2001; Gill 
et al. 2004b; Giocoli, Tormen & van den Bosch 2008). 

This phenomenon is further illustrated in Fig. [3] which 
shows the ratio between the final mass of a surviving satel- 
lite at the end of the simulation and its mass at the infall 
or orbit entry time as a function of infall time, i.e. the time 
it crosses the virial radius of the host for the first time. The 
tendency for mass loss is clearly seen as the ratios are typi- 
cally much lower than one especially for those haloes which 
became satellites early. Figure [4] shows this process in detail 
for a single example subhalo. During its evolution it first 
gains mass as the hierarchical formation scenario predicts, 
but then reaches a point at which it starts losing mass. As 
expected, the moment when the mass trend reverses is close 
to the time at which the subhalo becomes a satellite of a 
larger halo, which will be its host for the rest of the simula- 



tion. From that moment on the satellite is being stripped by 
tidal forces of the host. The mass loss obviously depends on 
the infall time but also on other parameters like the orbit. 
The presented example is a rather extreme case. According 
to Fig. [3] most haloes do not lose that much mass. Also not 
all subhaloes follow this path of evolution. Many of them 
just fall into the host halo and merge with it. 

Figure [5] shows the survival time of a satellite in terms 
of the number of pericentre passages it sustains. In the up- 
per panel only those satellites were selected which did not 
survive until the end of the simulation and were completely 
destroyed. Many subhaloes do not survive even a single peri- 
centre passage, they fall straight inside and merge with the 
host. This behaviour is in agreement with other A-body 
studies where it was found that majority of stars in the 
present MW halo comes from the most massive subhaloes 
that were accreted in the past (e.g. Bullock & Johnston 2005; 
Sales et al. 2007b). There are however still more than a hun- 
dred of satellites per halo, which were able to survive up to 
eight pericentre passages. Their orbital history is illustrated 
by the histogram in the lower panel of Fig. [5] Presumably, 
these might be the dwarf spheroidal galaxies of the LG that 
we presently see as, according to the tidal stirring scenario 
(Mayer et al. 2001; Klimentowski et al. 2009), a significant 
number of pericentric passages accompanied by substantial 
mass loss are required to transform a disky stellar compo- 
nent into a spheroid. 

Actually, at the end of the simulation, M31 has 198 and 
the MW 136 surviving subhaloes inside the virial radius. It 
is interesting to look at the distribution of the shapes of the 
orbits of these subhaloes. Figure[6]plots the apocentre versus 
pericentre (r a and r p ) for those subhaloes that survived until 
the end and completed at least one orbit around their host. 
Note that M31 subhaloes are strongly underrepresented be- 
cause many of them found inside the virial radius at the 
end have not yet completed a single orbit (see the next sec- 
tion). Figure [7] presents the distribution of the ratio r p /r-a 
of the same subhaloes. We can see that there is a strong 
preference for rather eccentric orbits in agreement with the 
results found by Gill et al. (2004b) and Diemand, Kuhlen & 
Madau (2007). 

Figure [8] shows the mass of subhaloes at the moment 
of orbit entry as a function of redshift at which the infall of 
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Figure 5. The survival of subhaloes. The upper panel shows 
the numbers of pericentre passages survived by subhaloes before 
their destruction therefore it shows those which eventually did not 
survive until the present. Subhaloes with zero pericentre passages 
were accreted directly on to the host halo and destroyed during 
their first passage. The lower panel corresponds to those subhaloes 
which survived until the end of the simulation and shows the 
fraction of subhaloes that survived a given number of pericentre 
passages to all surviving subhaloes separately in the M31 and 
in the MW. Here subhaloes with zero pericentre passages have 
entered their host recently and have not yet reached the pericentre 
of their orbit. In each panel the horizontal shading corresponds 
to the M31, the skewed one to the MW. 



the subhalo takes place. The Figure allows us to hypothe- 
size on the possible evolution of the subhaloes in different 
mass ranges. The majority of large mass subhaloes (with 
masses around 10 9 Mq) enter their orbits at around z = 1. 
These are most probably the progenitors of dwarf spheroidal 
galaxies. As shown in Klimentowski et. al (2009) a halo of 
this mass possessing a stellar disk at z = 1 has enough time 




pericenter [h _1 Mpc] 



Figure 6. The distribution of orbits for subhaloes that survived 
until the present, were found inside the virial radius of the host 
and completed at least one orbit around it. The filled and open 
circles show the apocentre versus pericentre distance respectively 
for M31 and MW subhaloes. The line indicates the circular orbits. 

20 I — i — i — i — | — i — i — i — | — i — i — i — | — i — i — i — | — i — i — i — 



15 



Z 10 



5 



0.2 0.4 0.6 0.8 1 

r p/ r a 

Figure 7. The distribution of peri- to apocentre ratio for the 
same sample of M31 and MW subhaloes as shown in Fig. [6] 

to form a dwarf spheroidal by the present by tidal stirring. 
The mass loss that occurs in this process can be very sig- 
nificant as demonstrated by Fig. [3] (cf. Mayer et al. 2001; 
Hayashi et al. 2003; Gill et al. 2004b; Kazantzidis, Moore & 
Mayer 2004b; Kampakoglou & Benson 2007; Klimentowski 
et al. 2007, 2009; Pefiarrubia, Navarro & McConnachie 2008; 
Giocoli et al. 2008) which at the end of the evolution leaves 
us with masses of the order of few times 1O 7 M0, as indeed 
measured for most dwarf spheroidals (e.g. Lokas et al. 2008; 
Lokas 2009). 

The subhaloes that start orbiting their host later most 
probably correspond to dwarf irregular galaxies, as they al- 
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Figure 8. The maximum mass of a subhalo as a function of red- 
shift at which this mass was reached. The mass corresponds to 
the moment when a subhalo enters its orbit around the host halo 
and starts losing mass, as shown in Fig. [4] The upper diagram is 
for the subhaloes of M31 while the lower one for the subhaloes of 
the MW. Only subhaloes which survived until the present time 
were plotted. The solid line in each panel indicates our complete- 
ness limit of 100 particles, corresponding to M = 2.5 X 10 7 h~ 1 
M . 



ready have the correct masses around 10 8 - 10 9 M Q and 
not enough time to evolve anymore. If they are also dwarf 
spheroidals today they would need to be dwarf spheroidals 
already at the moment of orbit entry. A process which could 
lead to the formation of a dwarf spheroidal (or a dwarf el- 
liptical) at large distance from a host galaxy is a collision of 
its progenitor with another similar object. 

Another class of satellites that can be distinguished in 
Fig. [8] are the galaxies which enter their orbits early and 
possess low masses of the order of a few times 10 7 M©. 
They would then evolve by tidal stirring, like normal dwarf 
spheroidals, leading to similar objects but of much smaller 
masses. These could be the progenitors of the ultra faint 



dwarfs with expected present masses of around 10 M© and 
smaller. As demonstrated in Klimentowski et al. (2009) the 
tidal stripping alone does not significantly change the mass- 
to-light ratios so if they indeed possess high M/L values, 
they must be due to the processes involving the evolution 
of the baryonic component, in particular the gas dynamics. 
Indeed, their very low initial masses could suppress star for- 
mation in the early stages (e.g. White & Rees 1978; Kauff- 
mann, White & Guiderdoni 1993; Mayer et al. 2007). Then 
they would host very few stars from the beginning which 
explains why they possess so few stars at the present time 
(Haiman, Rees & Loeb 1997). 

Obviously, the mass resolution of the simulation de- 
scribed here is not sufficient to study the tidal stripping 
of satellites in detail. Figures [3] and 2] should therefore be 
considered only as indicative that this process is indeed 
present. In order to make reliable predictions on its effec- 
tiveness one needs to follow the evolution of dwarf galaxies 
with much higher resolution, which usually means evolving 
a single dwarf in a fixed potential of the host, as done e.g. in 
Klimentowski et al. (2007, 2009). Our purpose here is rather 
to investigate other possible scenarios for the formation of 
early type dwarfs which involve neighbours. The following 
sections are therefore devoted to the evolution of satellites 
in groups and the merging of satellites. 



4 SATELLITES IN GROUPS 

It has been recently proposed that the infall of subhaloes on 
to a large halo might occur in groups. The claims have been 
based both on A r -body simulations (e.g. Li & Helmi 2008; 
Angulo et al. 2009) and observational data (e.g. D'Onghia 
& Lake 2008; Metz et al 2009a). The distribution of dwarf 
satellites of the MW also suggests that they did not in- 
fall from completely random directions and form a disk- 
like structure (Metz et al. 2007, 2009a). If such groups of 
infalling haloes were bound the probabilities of their close 
interactions and mergers before the infall would be much 
higher than for randomly distributed haloes. In this case 
the haloes inside such groups could be the progenitors of 
at least some dwarf spheroidal or elliptical galaxies forming 
them either by tidal interactions or collisions. In this section 
we study such a possibility in detail. 



4.1 Group finding algorithm 

Our group finding algorithm is a generalized version of the 
friends-of-friends (FOF) algorithm (Davis et al. 1985). As in 
the classical FOF procedure, we define a linking length. An 
additional parameter is the number of simulation outputs or 
a time period. Two haloes are considered to form a group if 
they are closer to each other than the linking length / for a 
given time period p. Each pair of haloes is checked. If one 
of the haloes which meet the criteria is already a member 
of a group, then the other one is added. When both haloes 
are already members of different groups, then the groups are 
merged. The algorithm reduces to the standard FOF when 
the period is set to a single snapshot. 

The algorithm is capable of finding many different kinds 
of groups. We have not used all its possibilities though. One 
of the most important features of the algorithm is the fact 
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Figure 9. The distribution of group size (measured as number 
of group members) in M31 (upper panel) and the MW (lower 
panel). The groups were identified by applying the group finding 
algorithm with parameters: I = 100 h~ 1 kpc, p = 1.5 Gyr, t± = 7 
Gyr, t 2 = 10.5 Gyr for M31 and / = 100 h' 1 kpc, p = 1.5 Gyr, 
ti = 0.2 Gyr, t 2 = 13.4 Gyr for the MW. 



that a group does not have to be strictly defined in time. For 
example, two haloes could form a group at one time, while 
one of these haloes could form a group with a third halo at 
a different time. Although we have two different groups of 
two haloes they are linked together to form one group by 
a common halo. If we do not want such linking we need to 
define a time period (with the starting moment ti and the 
ending moment t 2 where ti — t\ > p), at which we want 
the algorithm to actually link the haloes. In practice, we use 
long time periods to find the groups and study their history, 
and then reapply the algorithm with a short time period at 
the time of infall to study the group behaviour in its host 
halo. The set of four parameters l,p,ti,t 2 completely define 
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Figure 10. The cumulative mass functions of the largest groups 
in M31 (LGA, solid line) and the MW (LGM, dashed line). 




Figure 11. M31 together with its largest subhalo visible to the 
upper right of the main halo. The size of the image is 200 ft~ 1 kpc. 



the algorithm. Obviously, it can be run on different sub- 
sets of haloes. We have run the algorithm for the subhaloes 
of the two most massive haloes, the M31 and the MW. A 
subhalo was defined as before as a halo which in at least 
one simulation output was closer to the host than one virial 
radius. 



4.2 Results of the group finding for M31 

We tested different sets of parameters of the algorithm. We 
were looking for a set of parameters which would recover 
whole complete groups but would also reject haloes that 
were members of the group for a short time. Results are quite 
stable and changing the linking length and the period by less 
than 30 per cent affects the group membership very little. 
Finally, two different sets of parameters of the group finder 
were used. The first run for M31 was made with I — 100 h^ 1 
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Figure 12. Trajectories of the subhaloes in the Large Group of 
Andromeda (LGA) starting around 6 Gyr from the beginning of 
the simulation. The upper panel shows the positions of subhaloes 
at 9 Gyr, the lower one at present (13.4 Gyr). In both panels the 
biggest halo is indicated by the thicker line. The dashed circle 
marks the present virial radius of M31. A square at the end of 
the line indicates that the halo was destroyed before the end of 
the simulation, a cross shows the position of a halo which survived 
until the present. 
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Figure 13. The directions of the angular momentum vectors of 
LGA subhaloes with respect to M31. The crosses indicate the 
angular coordinates of the vectors (in the spherical coordinates 
system) at the time of infall, t = 9 Gyr (upper panel) and at the 
last simulation output t = 13.4 Gyr (lower panel). The coordinate 
system was chosen so that the angular momentum vector of M31 
points along the z axis and the orientation of x and y (different in 
each panel) so that the data points appear centrally in the plots. 



kpc, p = 1.5 Gyr, ti = 0.2 Gyr, ti = 13.4 Gyr, thus in this 
case we look for groups in almost the whole duration of the 
simulation. The linking length of 100 h~ x kpc is sufficient 
for finding large groups of subhaloes. 

As a result we get one large group of 116 haloes and 
about 40 small groups with up to 10 haloes. It turns out 
that the large group is indeed a bound structure, while the 
smaller groups are just accidentally linked. The large group 
was the first target of the analysis. The next step was to 
narrow down the criteria of the group finder. In the second 
run we used I = 100 h~ 1 kpc, p — 1.5 Gyr, t\ = 7 Gyr, 
ti = 10.5 Gyr. The linking length remains the same, but 
the period at which the linking is made was reduced to the 
time between 7 and 10.5 Gyr which roughly corresponds to 



the infall time of the large group on to M31. This way we 
were able to study the large group itself at a time close to 
its accretion and later. The second run thus can be treated 
as a filter applied on the first run allowing to reject haloes 
which deviated from the group before infall. At this point 
the results are even less influenced by the exact values of the 
parameters. 

4.2.1 The large group of M31 

The upper panel of Fig. [9] shows the histogram of the dis- 
tribution of the group size following from the second run of 
the group finding algorithm for M31. With the adopted pa- 
rameters the largest group consists of 30 haloes. From now 
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so that the trajectories are continuous between = and 360°. Different symbols correspond to snapshots at different times (in Gyr) 
indicated in the legend. 



on we call it the Large Group of Andromeda (LGA). Few 
tens of smaller groups are also found. The cumulative mass 
function of the LGA shown with the solid line in Fig. [TU] 
illustrates its nature. The group consists of one very large 
halo of mass M = 7.5 x 10 10 h' 1 M , well visible as the 
largest subhalo in the image of M31 shown in Fig. 111! and 
several tens of smaller haloes which are actually its satellites. 
Figure [12] shows the spatial distribution of the haloes in the 
group at two different times. We did not find any mergers 
inside such a group. 



4-2.2 Observational effects of a LGA-like group 

We now consider the present effects of an infall of a LGA-like 
group in the past. In this case we have a large halo which 
falls on to M31 along with the set of its satellites. In the 
simulation the large halo survives until the end, but eventu- 
ally it will merge with M31. We are more interested in the 
properties of the satellites. The question is whether we can 
detect such a past merger by studying only the properties 
of dwarfs at present. 

Figure [13] shows the distribution of the angular mo- 
menta vectors of the LGA subhaloes with respect to M31. 
Only the direction of the vector is shown in the spherical co- 
ordinate system. The orientation of the system was chosen 
so that the z axis points along the total angular momentum 
of M31. One can clearly see that during the infall (upper 
panel) angular momenta of the subhaloes are correlated as 
postulated by Li & Helmi (2008) who defined groups by the 
separation of angular momenta vectors. In our case how- 
ever only some of the haloes show this pattern and several 
others have their angular momentum vectors differing more 
than the value of only 10° used by Li & Helmi. This sug- 
gests that the method of group finding by comparing only 
angular momentum orientations may miss larger objects or 
other bound subhaloes. The lower panel of Fig. [T3] shows 
the distribution in the last simulation output, i.e. 4.4 Gyr 
later. This time period roughly corresponds to half of the 
group orbit around M31, as can be seen in Fig. 1121 During 
that time the angular momenta were modified and no clear 
pattern seems to be preserved. 

Figure [14] shows the LGA haloes orbits projected on 



to the sky. The observer was placed in the M31 centre and 
the coordinate system was again chosen so that the z axis 
is aligned with the angular momentum of M31. The haloes 
are infalling from one well-defined region in the sky but the 
spread increases with time. It reaches a maximum around 
12.4 Gyr when the haloes occupy almost half of the sky. 
Then the haloes from the main group seem to fall back on the 
group again, but some other are left in the opposite part of 
the sky. This example suggests that the alignment of angular 
momenta of infalling satellites is not well conserved, even 
though the group has not yet decayed (see also Libeskind et 
al. 2007). One should thus be very careful when trying to 
reproduce the histories of dwarf galaxies using their present 
proper motions. Metz et al. (2007) studied the so-called disk 
of satellites in the LG. They found that the distribution of 
MW and M31 satellites is not isotropic, and it rather forms 
a disk. It has been claimed that this kind of structure could 
be an effect of a group infall on to the MW halo (Li & Helmi 
2008). Based on our result we conclude that such a disk is 
probably not an effect of a group infall unless it happened 
very recently. 



4.2.3 Smaller groups of M31 

The upper panel of Fig. [5] shows that apart from the large 
group there is a number of smaller groups consisting of usu- 
ally two, sometimes several haloes. It is interesting to check 
whether haloes in such groups could strongly interact with 
each other. We have calculated the escape velocities from 
such groups and compared them to the true velocities of the 
haloes. It turns out that for the majority of the groups the 
subhaloes are only by chance falling on to the host halo from 
the same direction and are not bound to each other. Several 
of such groups of two haloes are marginally bound for very 
short periods of time. We studied in detail mass histories of 
such cases and concluded that no significant mass transfer 
is present. These haloes are not massive enough to influence 
their neighbourhood by tidal interactions and to stay bound 
for a longer period of time. 
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4.3 Results of the group finding for MW 

The same algorithm was applied to the MW subhaloes. The 
analysis shows that there are no groups similar to the LGA. 
The lower panel of Fig. [9] shows the distribution of group 
sizes obtained with I = 100 h" 1 kpc, p = 1.5 Gyr, t\ = 0.2 
Gyr, t2 = 13.4 Gyr. We get a few tens of smaller groups 
and a larger one consisting of 23 subhaloes, which we will 
call the Large Group of MW or LGM (with these param- 
eters LGA had 116 members). A closer inspection shows 
that this group, although of different nature than LGA, is 
an interesting object by itself. Figure \W\ presents the cumu- 
lative mass function of this group with the dashed line. The 
group consists of four larger subhaloes with masses around 
M = lO 9 /*," 1 M and some smaller ones. Thus, in contrast to 
LGA, this is a group of satellites of comparable size, rather 
than a single large halo with its subhaloes. Unfortunately 
the infall time of the group is very late, it starts at around 
12 Gyr from the beginning of the simulation and is not com- 
pleted by the end. This prevents us from studying the future 
evolution of such a group inside the MW halo. 

The group was formed around 4 Gyr from the beginning 
of the simulation. During most of the time it consists of 12 
haloes, while the rest is accreted late, when the group falls 
on to the MW. A careful analysis shows that there is no sig- 
nificant mass transfer between large members. The smaller 
subhaloes tend to lose mass but this process is much slower 
than for MW subhaloes and in total it does not amount to 
more than a few percent of the initial mass. The larger sub- 
haloes tend to gain mass both from smaller satellites and 
from other infalling objects, but this process again is very 
slow and insignificant. We conclude that the group members 
do not undergo any significant evolution inside the group. 



5 MERGERS AND INTERACTIONS OF THE 
HALOES 

In the previous section we have shown that being a mem- 
ber of a group of subhaloes does not significantly affect the 
evolution of a given subhalo and thus cannot by itself lead 
to the formation of an early type galaxy. Another possible 
channel by which such objects could form are interactions 
and direct mergers of subhaloes. This issues have already 
been addressed to some extent by Knebe, Gill & Gibson 
(2004) and Knebe et al. (2006) where they found that on 
average 30 per cent of the substructure population experi- 
enced encounters and that such interactions can account for 
a significant fraction of mass loss. In this section we address 
the question of how often such mergers and interactions oc- 
cur and whether they can indeed lead to the formation of 
dwarf spheroidals and ellipticals. 
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Figure 15. The distribution of interactions for the MW and M31 
surviving subhaloes in time for the parameters mi = 0.1, ma = 
0.3. The number of interactions was expressed as a fraction of all 
interactions and the histogram was normalized to unity. 



For example, the values mi = 0.5, rrt2 = 0.3 mean that a 
minimal interaction would be with haloes of mass equal to 
half the mass of the studied halo from which it gains 30 
per cent. For an interaction to be interesting, m\ has to be 
large enough, while m,2 needs to be some significant fraction 
of mass of the other haloes. It is worth noting that a halo 
might become a subhalo of a larger halo, but then leave it 
again (as found in cosmological simulations by e.g. Gill et al. 
2005; Warnick, Knebe & Power 2008; Ludlow et al. 2009). In 
such a case the algorithm is capable of finding the interac- 
tion for both haloes even though it is not a classical merger. 
The results will depend on the assumed parameters and the 
amount of matter exchanged. 

The motivation of this approach is the following. In 
most studies only major mergers are considered. However 
here we are interested in the probability for a halo to po- 
tentially form a dwarf spheroidal or elliptical galaxy by an 
interaction. Usually, a strong interaction between two haloes 
would eventually lead to a merger anyway. But the merger 
itself could happen between haloes which have already ex- 
changed large amounts of mass in past interactions. The 
resulting merger could seem not to be significant. Our ap- 
proach allows us to detect each of those interactions and 
find the strongest one instead of studying only the merger 
itself. We would like to know what was the magnitude of the 
strongest of these interactions, as this one had the largest 
chance to transform the galaxy. 



5.1 Algorithm for finding interacting haloes 

The algorithm presented here was based on finding and se- 
lecting interacting haloes. We look for haloes between which 
particles are exchanged. For each halo an interaction is de- 
fined by the following pair of parameters: the mass of all 
the other haloes taking part in the encounter expressed as a 
fraction of the studied halo's mass, mi, and the mass frac- 
tion gained by the studied halo during the encounter, m%. 



5.2 Interactions of the infalling haloes 

We consider only the interactions occurring for those haloes 
which are still on their way to the host and have not yet be- 
come satellites. In this analysis we include only those haloes 
of the MW and M31 which survived until the end of the 
simulation with 100 particles or more, as only those are mas- 
sive enough to be interesting from the observational point of 
view. We check all detected events so a given surviving sub- 
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Figure 16. The dependence of the number of interactions on 
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that happened to the surviving sample of 198 M31 and 136 MW 
subhaloes as a function of m\ . The different line types correspond 
to different values of mi as indicated by the legend. The numbers 
give the ratio of the number of interactions to the total number 
of subhaloes. 



halo could have had more that one interaction in the past. 
We also require that at the moment of interaction the stud- 
ied halo needs to have at least 50 particles to avoid numerical 
noise. 

As already mentioned, at the end of the simulation 
M31 has 198 and the MW 136 surviving subhaloes of mass 
M > 2.5 x 10 7 h~ x Mq. For example, if we assume that the 
set of parameters mi = 0.5, mi — 0.3 would be characteris- 
tic for a major interaction and apply it, then the algorithm 
detects a total number of 28 interactions for both M31 and 
MW subhaloes. This means that less than 10 percent of sur- 
viving satellites had a major interaction in the past. This is 
a significant number, but not high enough to explain the ori- 
gin of all dwarf spheroidal or elliptical satellites. Figure [15] 
shows a combined distribution of both M31 and MW inter- 
actions in time. Most of the mergers occurred quite early, 
around redshift z = 2 corresponding to 3 Gyr since begin- 
ning of the simulation. Such a collision would not have to 
produce an early type dwarf galaxy immediately, since it 
would still have a lot of time to evolve in the tidal field of 
its host. 

Note that the shape of the distribution shown in Fig. [15] 
does not depend strongly on the assumed parameters mi 
and mi; the interactions still occur rather early if we vary 
those parameters. The total number of interactions detected 
depends however on these parameters. Figure [16] presents 
the results of exploring the parameter space in terms of the 
fraction of subhaloes that experienced interaction (or a num- 
ber of interactions per subhalo) with a given mi and mi. 
As expected, minor interactions are more numerous and the 
fraction of affected subhaloes decreases strongly with grow- 
ing mi and mi. 

It turns out that a significant fraction (15-25 percent, 



taking the smallest mi and the whole range of mi in Fig. 1161) 
of subhaloes had some kind of interaction in the past, but 
only several percent had a strong interaction that could be 
called a major merger. This suggests that the mergers, while 
present in the histories of some dwarf galaxies, could not 
explain the large numbers of spheroidals and ellipticals in 
the LG. Still, majority of dwarf galaxies did not have even 
a slight interaction with another dwarf in the past. 

However, this result depends on the assumption that we 
count the interactions only when the halo has at least 50 par- 
ticles, corresponding to the mass of M = 1.25 x lO 7 /^ 1 M Q , 
which we find to be the lowest limit at which this study still 
makes sense. We have verified that the assumed threshold 
of 50 particles for the algorithm is rather realistic: with in- 
creasing threshold the numbers of interactions decrease very 
slowly, while lowering the threshold would produce many 
more interactions. This could be easily understood: accord- 
ing to Fig. ll5l most of the encounters are expected to happen 
early when the haloes are in general less massive. Interac- 
tions with low particle numbers may however be affected by 
numerical noise and therefore their statistics is not reliable. 



6 DISCUSSION 

We have analyzed a constrained simulation of the local Uni- 
verse which reproduces the main properties of the LG. The 
simulation outputs were used to study the population of 
subhaloes around the largest galaxies. Our attention was 
focused on those signatures of the evolution that can shed 
some light on the possible scenarios of the formation of early 
type dwarf galaxies, dwarf ellipticals and dwarf spheroidals, 
that we presently observe. Assuming that the progenitors of 
these objects are baryonic disks embedded in dark matter 
haloes, some mechanisms are required that could transform 
the disks into bulges. In the case of purely gravitational in- 
teractions, such an evolution can occur by three channels: 
the tidal stirring in the gravitational field of the host galaxy, 
interaction of small haloes infalling together in groups and 
mergers between prospective subhaloes. 

A quantitative study of the effect of tidal forces on a 
single subhalo usually requires a different simulation setup 
and much higher resolution (e.g. Klimentowski et al. 2007, 
2009). In the present context the tidal stripping manifests 
itself in the statistical properties of our subhalo population. 
We found that the mass function of satellites evolves signifi- 
cantly between z = 1 and z — so that the number of small 
subhaloes increases and the number of larger ones decreases, 
as expected if subhaloes lose mass by tides. On average, the 
mass of satellites at the end of the simulation is smaller than 
at the moment of their infall on to the halo of their host. For 
most subhaloes the transition from the phase of mass accre- 
tion to the phase of mass loss can be easily identified and 
corresponds to the moment of entry into the host halo. We 
find that this maximum mass or entrance moment can occur 
at a wide range of redshifts but typically more massive sub- 
haloes enter the vicinity of their hosts earlier, around z = 1 
and before (see also Diemand et al. 2007). This distribution 
in redshift suggests that those massive subhaloes could be 
the progenitors of present-day dwarf spheroidals since they 
have still a lot of time to evolve. The less massive subhaloes 
entering later could still be dwarf irregulars today on one 
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of their first passages around the MW like the Magellanic 
Clouds (Besla et al. 2007). 

We also studied the evolution of satellites in groups. Al- 
though many groups were identified both around M31 and 
the MW, they are usually loosely bound and their member 
subhaloes do not strongly interact with each other. There- 
fore we conclude that being a member of a group cannot re- 
sult in any morphological transformation of a dwarf galaxy. 
Membership in a group could however explain some par- 
ticular distributions of satellites around big galaxies, pro- 
vided that the group infall occurred recently. Using a large 
group identified in the vicinity of our simulated M31 we have 
demonstrated that the group forms a coherent structure at 
infall but dissipates on a rather short timescale. 

Mergers between subhaloes offer another possibility for 
the formation of early type dwarfs. We found that around 10 
percent of present-day subhaloes of M31 and the MW un- 
derwent a major interaction with another halo in the past. 
Most of the events took place early on, around 7—11 Gyr 
ago, when the objects have not yet become satellites of big 
galaxies. This suggests that dwarf ellipticals could be a prod- 
uct of mergers rather than tidal evolution, as they tend to 
be more isolated objects found further from the MW than 
dwarf spheroidals. 
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